Libraries and paths
# Libraries
library(easypackages)
libraries("here","tidyverse","ggeasy","patchwork","reshape2")
fontSize = 16
dotSize = 2
alphaVal = 0.5
dot_cols2use = c("#404040","#a3142e")
dtype_cols2use = c("#7ab5ed","#171c9e")
# Paths
codedir = here("code")
plotdir = here("plots")
datadir = here("data")
# functions that will be used over and over again
# make_scatterplot
make_scatterplot <- function(data2plot, x_var, y_var, strat_var,
xlabel2use, ylabel2use, title2use,
vlineval=NULL,plot_fit_line=FALSE,plot_diag_line=FALSE,
dotSize=2, fontSize=16, alphaVal=0.5,
cols2use=c("#404040","#a3142e")
){
p = ggplot(data = data2plot, aes_string(x = x_var, y = y_var, colour= strat_var))
if (is.null(alphaVal)){
p = p + geom_point(size = dotSize)
} else{
p = p + geom_point(size = dotSize, alpha = alphaVal)
}
p = p +
scale_colour_manual(values=cols2use) +
xlab(xlabel2use) +
ylab(ylabel2use) +
ggtitle(title2use) +
easy_center_title() +
guides(colour="none")+
theme_bw() +
theme(plot.title = element_text(hjust = 0.5, size=fontSize),
strip.text.y = element_text(size=fontSize),
axis.text.x = element_text(size=fontSize),
axis.text.y = element_text(size=fontSize-2),
axis.title.x = element_text(size=fontSize),
axis.title.y = element_text(size=fontSize),
strip.text = element_text(size = fontSize))
if (!is.null(vlineval)){
p = p + geom_vline(xintercept = vlineval)
}
if (plot_fit_line){
p = p + geom_smooth(method = lm, se=FALSE, aes_string(group=strat_var))
}
if (plot_diag_line){
p = p + geom_abline(intercept = 0, slope = 1, colour = "black")
}
return(p)
} # function make_scatterplot
# out_of_sample_prediction
out_of_sample_prediction <- function(trn_data, val_data, dv_var, pred_var){
form2use = as.formula(sprintf("%s ~ %s", dv_var, pred_var))
mod2use = lm(formula = form2use, data = trn_data)
val_data$pred = predict(mod2use, val_data)
# compute R-squared
rss = sum((val_data[,"pred"] - val_data[,dv_var])^2) ## residual sum of squares
tss = sum((val_data[,dv_var] - mean(val_data[,dv_var]))^2) ## total sum of squares
rsq = 1 - (rss/tss)
res_list = list(trn_data = trn_data,
val_data = val_data,
mod2use = mod2use,
rsq = rsq)
return(res_list)
} # out_of_sample_prediction
Read in data
g_eeg_H_data = read.csv(file.path(datadir,"tidy_H_g_EEG_data.csv"))
g_lfp_H_data = read.csv(file.path(datadir,"tidy_H_g_LFP_data.csv"))
g_h_data = rbind(g_eeg_H_data, g_lfp_H_data)
g_eeg_slope_data = read.csv(file.path(datadir,"tidy_slope_g_EEG_data.csv"))
g_lfp_slope_data = read.csv(file.path(datadir,"tidy_slope_g_LFP_data.csv"))
g_slope_data = rbind(g_eeg_slope_data, g_lfp_slope_data)
g_eeg_pow_data = read.csv(file.path(datadir,"tidy_pow_g_EEG_data.csv"))
g_lfp_pow_data = read.csv(file.path(datadir,"tidy_pow_g_LFP_data.csv"))
g_pow_data = rbind(g_eeg_pow_data, g_lfp_pow_data)
camk_eeg_H_data = read.csv(file.path(datadir,"tidy_H_camk_EEG_data.csv"))
camk_lfp_H_data = read.csv(file.path(datadir,"tidy_H_camk_LFP_data.csv"))
camk_h_data = rbind(camk_eeg_H_data, camk_lfp_H_data)
camk_eeg_slope_data = read.csv(file.path(datadir,"tidy_slope_camk_EEG_data.csv"))
camk_lfp_slope_data = read.csv(file.path(datadir,"tidy_slope_camk_LFP_data.csv"))
camk_slope_data = rbind(camk_eeg_slope_data, camk_lfp_slope_data)
camk_eeg_pow_data = read.csv(file.path(datadir,"tidy_pow_camk_EEG_data.csv"))
camk_lfp_pow_data = read.csv(file.path(datadir,"tidy_pow_camk_LFP_data.csv"))
camk_pow_data = rbind(camk_eeg_pow_data, camk_lfp_pow_data)
hm4di_eeg_H_data = read.csv(file.path(datadir,"tidy_H_hm4di_EEG_data.csv"))
hm4di_lfp_H_data = read.csv(file.path(datadir,"tidy_H_hm4di_LFP_data.csv"))
hm4di_h_data = rbind(hm4di_eeg_H_data, hm4di_lfp_H_data)
hm4di_eeg_slope_data = read.csv(file.path(datadir,"tidy_slope_hm4di_EEG_data.csv"))
hm4di_lfp_slope_data = read.csv(file.path(datadir,"tidy_slope_hm4di_LFP_data.csv"))
hm4di_slope_data = rbind(hm4di_eeg_slope_data, hm4di_lfp_slope_data)
hm4di_eeg_pow_data = read.csv(file.path(datadir,"tidy_pow_hm4di_EEG_data.csv"))
hm4di_lfp_pow_data = read.csv(file.path(datadir,"tidy_pow_hm4di_LFP_data.csv"))
hm4di_pow_data = rbind(hm4di_eeg_pow_data, hm4di_lfp_pow_data)
Plot Firing Rate by g
# Scatterplot of Firing rate by g
data2plot = g_h_data %>%
filter(dset=="training") %>%
filter(measure=="EEG") %>%
select(fr_exc, fr_inh, param) %>%
distinct() %>%
melt(id.vars=c("param"))
vlineval = 0.0885
yLabel = "Firing Rate (Hz)"
xLabel = expression("g =" ~ g[E]/g[I])
title2use = xLabel
cols2use = dot_cols2use
FR_g_plot = make_scatterplot(data2plot = data2plot,
x_var = "param",
y_var = "value",
strat_var = "variable",
vlineval = vlineval,
xlabel2use = xLabel,
ylabel2use = yLabel,
title2use = title2use,
cols2use = cols2use,
alphaVal = NULL)
FR_g_plot
# correlation with g in E neurons
data2use = data2plot %>% filter(variable=="fr_exc")
res = cor.test(data2use[,"param"], data2use[,"value"])
res
##
## Pearson's product-moment correlation
##
## data: data2use[, "param"] and data2use[, "value"]
## t = 19.051, df = 68, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.8704428 0.9482165
## sample estimates:
## cor
## 0.9177152
res$p.value
## [1] 5.698427e-29
# correlation with g in I neurons
data2use = data2plot %>% filter(variable=="fr_inh")
res = cor.test(data2use[,"param"], data2use[,"value"])
res
##
## Pearson's product-moment correlation
##
## data: data2use[, "param"] and data2use[, "value"]
## t = 67.452, df = 68, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.9880972 0.9954156
## sample estimates:
## cor
## 0.9926099
res$p.value
## [1] 5.032159e-64
Plot Firing Rate by Camk El
# Scatterplot of Firing rate by Camk El
data2plot = camk_h_data %>%
filter(dset=="training") %>%
filter(measure=="EEG") %>%
select(fr_exc, fr_inh, param) %>%
distinct() %>%
melt(id.vars=c("param"))
vlineval = -70
yLabel = "Firing Rate (Hz)"
xLabel = "Resting Potential (mV)"
title2use = "CamkII-hM3D(Gq)"
cols2use = dot_cols2use
FR_camk_El_plot = make_scatterplot(data2plot = data2plot,
x_var = "param",
y_var = "value",
strat_var = "variable",
vlineval = vlineval,
xlabel2use = xLabel,
ylabel2use = yLabel,
title2use = title2use,
cols2use = cols2use,
alphaVal = NULL)
FR_camk_El_plot
# correlation with g in E neurons
data2use = data2plot %>% filter(variable=="fr_exc")
res = cor.test(data2use[,"param"], data2use[,"value"])
res
##
## Pearson's product-moment correlation
##
## data: data2use[, "param"] and data2use[, "value"]
## t = 52.206, df = 49, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.984434 0.994953
## sample estimates:
## cor
## 0.9911301
res$p.value
## [1] 1.289149e-44
# correlation with g in I neurons
data2use = data2plot %>% filter(variable=="fr_inh")
res = cor.test(data2use[,"param"], data2use[,"value"])
res
##
## Pearson's product-moment correlation
##
## data: data2use[, "param"] and data2use[, "value"]
## t = 59.298, df = 49, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.9878897 0.9960781
## sample estimates:
## cor
## 0.9931044
res$p.value
## [1] 2.760435e-47
Plot Firing Rate by hM4Di El
# Scatterplot of Firing rate by hM4Di El
data2plot = hm4di_h_data %>%
filter(dset=="training") %>%
filter(measure=="EEG") %>%
select(fr_exc, fr_inh, param) %>%
distinct() %>%
melt(id.vars=c("param"))
vlineval = -70
yLabel = "Firing Rate (Hz)"
xLabel = "Resting Potential (mV)"
title2use = "hM4Di"
cols2use = dot_cols2use
FR_hm4di_El_plot = make_scatterplot(data2plot = data2plot,
x_var = "param",
y_var = "value",
strat_var = "variable",
vlineval = vlineval,
xlabel2use = xLabel,
ylabel2use = yLabel,
title2use = title2use,
cols2use = cols2use,
alphaVal = NULL)
FR_hm4di_El_plot
# correlation with g in E neurons
data2use = data2plot %>% filter(variable=="fr_exc")
res = cor.test(data2use[,"param"], data2use[,"value"])
res
##
## Pearson's product-moment correlation
##
## data: data2use[, "param"] and data2use[, "value"]
## t = 5.5862, df = 9, p-value = 0.0003402
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.5962093 0.9688472
## sample estimates:
## cor
## 0.8809954
res$p.value
## [1] 0.0003402203
# correlation with g in I neurons
data2use = data2plot %>% filter(variable=="fr_inh")
res = cor.test(data2use[,"param"], data2use[,"value"])
res
##
## Pearson's product-moment correlation
##
## data: data2use[, "param"] and data2use[, "value"]
## t = 17.813, df = 9, p-value = 2.515e-08
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.9456054 0.9965087
## sample estimates:
## cor
## 0.9861128
res$p.value
## [1] 2.5149e-08
Analysis of E:I ratio simulation
# ==============================================================================
# H
# Scatterplot against g
data2plot = g_h_data %>% filter(dset=="training")
vlineval = 0.0885
yLabel = "H"
xLabel = expression("g =" ~ g[E]/g[I])
title2use = xLabel
cols2use = dtype_cols2use
plot_fit_line = TRUE
H_g_plot = make_scatterplot(data2plot = data2plot,
x_var = "param",
y_var = "value",
strat_var = "measure",
vlineval = vlineval,
xlabel2use = xLabel,
ylabel2use = yLabel,
title2use = title2use,
cols2use = cols2use,
plot_fit_line = plot_fit_line)
H_g_plot
# correlation with g in LFP
data2use = data2plot %>% filter(measure=="LFP")
res = cor.test(data2use[,"param"], data2use[,"value"])
res
##
## Pearson's product-moment correlation
##
## data: data2use[, "param"] and data2use[, "value"]
## t = -11.035, df = 68, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## -0.8719076 -0.6973305
## sample estimates:
## cor
## -0.801049
res$p.value
## [1] 8.36119e-17
# correlation with g in EEG
data2use = data2plot %>% filter(measure=="EEG")
res = cor.test(data2use[,"param"], data2use[,"value"])
res
##
## Pearson's product-moment correlation
##
## data: data2use[, "param"] and data2use[, "value"]
## t = -11.365, df = 68, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## -0.8774729 -0.7093028
## sample estimates:
## cor
## -0.8093794
res$p.value
## [1] 2.262083e-17
# ==============================================================================
# out of sample prediction
trn_eeg_data = g_h_data %>% filter(dset=="training") %>% filter(measure=="EEG")
val_eeg_data = g_h_data %>% filter(dset=="validation") %>% filter(measure=="EEG")
oos_res_eeg = out_of_sample_prediction(trn_data = trn_eeg_data,
val_data = val_eeg_data,
dv_var="fr_tot",
pred_var="value")
print(sprintf("EEG R-sq = %f",oos_res_eeg$rsq))
## [1] "EEG R-sq = 0.558346"
trn_lfp_data = g_h_data %>% filter(dset=="training") %>% filter(measure=="LFP")
val_lfp_data = g_h_data %>% filter(dset=="validation") %>% filter(measure=="LFP")
oos_res_lfp = out_of_sample_prediction(trn_data = trn_lfp_data,
val_data = val_lfp_data,
dv_var="fr_tot",
pred_var="value")
print(sprintf("LFP R-sq = %f",oos_res_lfp$rsq))
## [1] "LFP R-sq = 0.671406"
# make plot of actual vs predicted
data2plot = rbind(oos_res_eeg$val_data, oos_res_lfp$val_data)
vlineval = NULL
yLabel = "Predicted Firing Rate (Hz)"
xLabel = "Actual Firing Rate (Hz)"
title2use = expression("g =" ~ g[E]/g[I])
cols2use = dtype_cols2use
plot_diag_line = TRUE
H_g_oos_plot = make_scatterplot(data2plot = data2plot,
x_var = "fr_tot",
y_var = "pred",
strat_var = "measure",
vlineval = vlineval,
xlabel2use = xLabel,
ylabel2use = yLabel,
title2use = title2use,
cols2use = cols2use,
plot_diag_line = plot_diag_line)
H_g_oos_plot
# ==============================================================================
# 1/f slope
# Scatterplot against g
data2plot = g_slope_data %>% filter(dset=="training")
vlineval = 0.0885
yLabel = "1/f slope"
xLabel = expression("g =" ~ g[E]/g[I])
title2use = xLabel
cols2use = dtype_cols2use
plot_fit_line = TRUE
slope_g_plot = make_scatterplot(data2plot = data2plot,
x_var = "param",
y_var = "value",
strat_var = "measure",
vlineval = vlineval,
xlabel2use = xLabel,
ylabel2use = yLabel,
title2use = title2use,
cols2use = cols2use,
plot_fit_line = plot_fit_line)
slope_g_plot
# correlation with g in LFP
data2use = data2plot %>% filter(measure=="LFP")
res = cor.test(data2use[,"param"], data2use[,"value"])
res
##
## Pearson's product-moment correlation
##
## data: data2use[, "param"] and data2use[, "value"]
## t = 16.573, df = 68, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.8362499 0.9338236
## sample estimates:
## cor
## 0.8953001
res$p.value
## [1] 1.411939e-25
# correlation with g in EEG
data2use = data2plot %>% filter(measure=="EEG")
res = cor.test(data2use[,"param"], data2use[,"value"])
res
##
## Pearson's product-moment correlation
##
## data: data2use[, "param"] and data2use[, "value"]
## t = 17.487, df = 68, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.8501953 0.9397321
## sample estimates:
## cor
## 0.9044783
res$p.value
## [1] 7.280804e-27
# ==============================================================================
# out of sample prediction
trn_eeg_data = g_slope_data %>% filter(dset=="training") %>% filter(measure=="EEG")
val_eeg_data = g_slope_data %>% filter(dset=="validation") %>% filter(measure=="EEG")
oos_res_eeg = out_of_sample_prediction(trn_data = trn_eeg_data,
val_data = val_eeg_data,
dv_var="fr_tot",
pred_var="value")
print(sprintf("EEG R-sq = %f",oos_res_eeg$rsq))
## [1] "EEG R-sq = 0.747769"
trn_lfp_data = g_slope_data %>% filter(dset=="training") %>% filter(measure=="LFP")
val_lfp_data = g_slope_data %>% filter(dset=="validation") %>% filter(measure=="LFP")
oos_res_lfp = out_of_sample_prediction(trn_data = trn_lfp_data,
val_data = val_lfp_data,
dv_var="fr_tot",
pred_var="value")
print(sprintf("LFP R-sq = %f",oos_res_lfp$rsq))
## [1] "LFP R-sq = 0.711435"
# make plot of actual vs predicted
data2plot = rbind(oos_res_eeg$val_data, oos_res_lfp$val_data)
vlineval = NULL
yLabel = "Predicted Firing Rate (Hz)"
xLabel = "Actual Firing Rate (Hz)"
title2use = expression("g =" ~ g[E]/g[I])
cols2use = dtype_cols2use
plot_diag_line = TRUE
slope_g_oos_plot = make_scatterplot(data2plot = data2plot,
x_var = "fr_tot",
y_var = "pred",
strat_var = "measure",
vlineval = vlineval,
xlabel2use = xLabel,
ylabel2use = yLabel,
title2use = title2use,
cols2use = cols2use,
plot_diag_line = plot_diag_line)
slope_g_oos_plot
# ==============================================================================
# Power
# Scatterplot against g
data2plot = g_pow_data %>% filter(dset=="training")
vlineval = 0.0885
yLabel = "Broadband Power (a.u.)"
xLabel = expression("g =" ~ g[E]/g[I])
title2use = xLabel
cols2use = dtype_cols2use
plot_fit_line = TRUE
pow_g_plot = make_scatterplot(data2plot = data2plot,
x_var = "param",
y_var = "value",
strat_var = "measure",
vlineval = vlineval,
xlabel2use = xLabel,
ylabel2use = yLabel,
title2use = title2use,
cols2use = cols2use,
plot_fit_line = plot_fit_line)
pow_g_plot
# correlation with g in LFP
data2use = data2plot %>% filter(measure=="LFP")
res = cor.test(data2use[,"param"], data2use[,"value"])
res
##
## Pearson's product-moment correlation
##
## data: data2use[, "param"] and data2use[, "value"]
## t = 34.722, df = 68, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.9566735 0.9831489
## sample estimates:
## cor
## 0.9729376
res$p.value
## [1] 5.370013e-45
# correlation with g in EEG
data2use = data2plot %>% filter(measure=="EEG")
res = cor.test(data2use[,"param"], data2use[,"value"])
res
##
## Pearson's product-moment correlation
##
## data: data2use[, "param"] and data2use[, "value"]
## t = 34.216, df = 68, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.9554483 0.9826658
## sample estimates:
## cor
## 0.9721658
res$p.value
## [1] 1.37964e-44
# ==============================================================================
# out of sample prediction
trn_eeg_data = g_pow_data %>% filter(dset=="training") %>% filter(measure=="EEG")
val_eeg_data = g_pow_data %>% filter(dset=="validation") %>% filter(measure=="EEG")
oos_res_eeg = out_of_sample_prediction(trn_data = trn_eeg_data,
val_data = val_eeg_data,
dv_var="fr_tot",
pred_var="value")
print(sprintf("EEG R-sq = %f",oos_res_eeg$rsq))
## [1] "EEG R-sq = 0.911498"
trn_lfp_data = g_pow_data %>% filter(dset=="training") %>% filter(measure=="LFP")
val_lfp_data = g_pow_data %>% filter(dset=="validation") %>% filter(measure=="LFP")
oos_res_lfp = out_of_sample_prediction(trn_data = trn_lfp_data,
val_data = val_lfp_data,
dv_var="fr_tot",
pred_var="value")
print(sprintf("LFP R-sq = %f",oos_res_lfp$rsq))
## [1] "LFP R-sq = 0.912992"
# make plot of actual vs predicted
data2plot = rbind(oos_res_eeg$val_data, oos_res_lfp$val_data)
vlineval = NULL
yLabel = "Predicted Firing Rate (Hz)"
xLabel = "Actual Firing Rate (Hz)"
title2use = expression("g =" ~ g[E]/g[I])
cols2use = dtype_cols2use
plot_diag_line = TRUE
pow_g_oos_plot = make_scatterplot(data2plot = data2plot,
x_var = "fr_tot",
y_var = "pred",
strat_var = "measure",
vlineval = vlineval,
xlabel2use = xLabel,
ylabel2use = yLabel,
title2use = title2use,
cols2use = cols2use,
plot_diag_line = plot_diag_line)
pow_g_oos_plot
Analysis of CamK
title2use = "CamkII-hM3D(Gq)"
# ==============================================================================
# H
# Scatterplot against El
data2plot = camk_h_data %>% filter(dset=="training")
vlineval = -70
yLabel = "H"
xLabel = "Resting Potential (mV)"
cols2use = dtype_cols2use
plot_fit_line = TRUE
H_camk_El_plot = make_scatterplot(data2plot = data2plot,
x_var = "param",
y_var = "value",
strat_var = "measure",
vlineval = vlineval,
xlabel2use = xLabel,
ylabel2use = yLabel,
title2use = title2use,
cols2use = cols2use,
plot_fit_line = plot_fit_line)
H_camk_El_plot
# correlation with El in LFP
data2use = data2plot %>% filter(measure=="LFP")
res = cor.test(data2use[,"param"], data2use[,"value"])
res
##
## Pearson's product-moment correlation
##
## data: data2use[, "param"] and data2use[, "value"]
## t = -37.543, df = 49, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## -0.9903432 -0.9703590
## sample estimates:
## cor
## -0.9830581
res$p.value
## [1] 9.037516e-38
# correlation with El in EEG
data2use = data2plot %>% filter(measure=="EEG")
res = cor.test(data2use[,"param"], data2use[,"value"])
res
##
## Pearson's product-moment correlation
##
## data: data2use[, "param"] and data2use[, "value"]
## t = -47.763, df = 49, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## -0.9939838 -0.9814634
## sample estimates:
## cor
## -0.9894305
res$p.value
## [1] 9.272795e-43
# ==============================================================================
# out of sample prediction
trn_eeg_data = camk_h_data %>% filter(dset=="training") %>% filter(measure=="EEG")
val_eeg_data = camk_h_data %>% filter(dset=="validation") %>% filter(measure=="EEG")
oos_res_eeg = out_of_sample_prediction(trn_data = trn_eeg_data,
val_data = val_eeg_data,
dv_var="fr_tot",
pred_var="value")
print(sprintf("EEG R-sq = %f",oos_res_eeg$rsq))
## [1] "EEG R-sq = 0.953605"
trn_lfp_data = camk_h_data %>% filter(dset=="training") %>% filter(measure=="LFP")
val_lfp_data = camk_h_data %>% filter(dset=="validation") %>% filter(measure=="LFP")
oos_res_lfp = out_of_sample_prediction(trn_data = trn_lfp_data,
val_data = val_lfp_data,
dv_var="fr_tot",
pred_var="value")
print(sprintf("LFP R-sq = %f",oos_res_lfp$rsq))
## [1] "LFP R-sq = 0.952495"
# make plot of actual vs predicted
data2plot = rbind(oos_res_eeg$val_data, oos_res_lfp$val_data)
vlineval = NULL
yLabel = "Predicted Firing Rate (Hz)"
xLabel = "Actual Firing Rate (Hz)"
cols2use = dtype_cols2use
plot_diag_line = TRUE
H_camk_El_oos_plot = make_scatterplot(data2plot = data2plot,
x_var = "fr_tot",
y_var = "pred",
strat_var = "measure",
vlineval = vlineval,
xlabel2use = xLabel,
ylabel2use = yLabel,
title2use = title2use,
cols2use = cols2use,
plot_diag_line = plot_diag_line)
H_camk_El_oos_plot
# ==============================================================================
# 1/f slope
# Scatterplot against El
data2plot = camk_slope_data %>% filter(dset=="training")
vlineval = -70
yLabel = "1/f slope"
xLabel = "Resting Potential (mV)"
cols2use = dtype_cols2use
plot_fit_line = TRUE
slope_camk_El_plot = make_scatterplot(data2plot = data2plot,
x_var = "param",
y_var = "value",
strat_var = "measure",
vlineval = vlineval,
xlabel2use = xLabel,
ylabel2use = yLabel,
title2use = title2use,
cols2use = cols2use,
plot_fit_line = plot_fit_line)
slope_camk_El_plot
# correlation with g in LFP
data2use = data2plot %>% filter(measure=="LFP")
res = cor.test(data2use[,"param"], data2use[,"value"])
res
##
## Pearson's product-moment correlation
##
## data: data2use[, "param"] and data2use[, "value"]
## t = 12.396, df = 49, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.7830842 0.9244918
## sample estimates:
## cor
## 0.8707545
res$p.value
## [1] 1.014247e-16
# correlation with g in EEG
data2use = data2plot %>% filter(measure=="EEG")
res = cor.test(data2use[,"param"], data2use[,"value"])
res
##
## Pearson's product-moment correlation
##
## data: data2use[, "param"] and data2use[, "value"]
## t = 13.587, df = 49, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.8123856 0.9353838
## sample estimates:
## cor
## 0.8889508
res$p.value
## [1] 3.062337e-18
# ==============================================================================
# out of sample prediction
trn_eeg_data = camk_slope_data %>% filter(dset=="training") %>% filter(measure=="EEG")
val_eeg_data = camk_slope_data %>% filter(dset=="validation") %>% filter(measure=="EEG")
oos_res_eeg = out_of_sample_prediction(trn_data = trn_eeg_data,
val_data = val_eeg_data,
dv_var="fr_tot",
pred_var="value")
print(sprintf("EEG R-sq = %f",oos_res_eeg$rsq))
## [1] "EEG R-sq = 0.864963"
trn_lfp_data = camk_slope_data %>% filter(dset=="training") %>% filter(measure=="LFP")
val_lfp_data = camk_slope_data %>% filter(dset=="validation") %>% filter(measure=="LFP")
oos_res_lfp = out_of_sample_prediction(trn_data = trn_lfp_data,
val_data = val_lfp_data,
dv_var="fr_tot",
pred_var="value")
print(sprintf("LFP R-sq = %f",oos_res_lfp$rsq))
## [1] "LFP R-sq = 0.849740"
# make plot of actual vs predicted
data2plot = rbind(oos_res_eeg$val_data, oos_res_lfp$val_data)
vlineval = NULL
yLabel = "Predicted Firing Rate (Hz)"
xLabel = "Actual Firing Rate (Hz)"
cols2use = dtype_cols2use
plot_diag_line = TRUE
slope_camk_El_oos_plot = make_scatterplot(data2plot = data2plot,
x_var = "fr_tot",
y_var = "pred",
strat_var = "measure",
vlineval = vlineval,
xlabel2use = xLabel,
ylabel2use = yLabel,
title2use = title2use,
cols2use = cols2use,
plot_diag_line = plot_diag_line)
slope_camk_El_oos_plot
# ==============================================================================
# Power
# Scatterplot against g
data2plot = camk_pow_data %>% filter(dset=="training")
vlineval = -70
yLabel = "Broadband Power (a.u.)"
xLabel = "Resting Potential (mV)"
cols2use = dtype_cols2use
plot_fit_line = TRUE
pow_camk_El_plot = make_scatterplot(data2plot = data2plot,
x_var = "param",
y_var = "value",
strat_var = "measure",
vlineval = vlineval,
xlabel2use = xLabel,
ylabel2use = yLabel,
title2use = title2use,
cols2use = cols2use,
plot_fit_line = plot_fit_line)
pow_camk_El_plot
# correlation with g in LFP
data2use = data2plot %>% filter(measure=="LFP")
res = cor.test(data2use[,"param"], data2use[,"value"])
res
##
## Pearson's product-moment correlation
##
## data: data2use[, "param"] and data2use[, "value"]
## t = 50.2, df = 49, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.9831882 0.9945468
## sample estimates:
## cor
## 0.9904176
res$p.value
## [1] 8.490014e-44
# correlation with g in EEG
data2use = data2plot %>% filter(measure=="EEG")
res = cor.test(data2use[,"param"], data2use[,"value"])
res
##
## Pearson's product-moment correlation
##
## data: data2use[, "param"] and data2use[, "value"]
## t = 48.374, df = 49, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.9819199 0.9941328
## sample estimates:
## cor
## 0.9896919
res$p.value
## [1] 5.036442e-43
# ==============================================================================
# out of sample prediction
trn_eeg_data = camk_pow_data %>% filter(dset=="training") %>% filter(measure=="EEG")
val_eeg_data = camk_pow_data %>% filter(dset=="validation") %>% filter(measure=="EEG")
oos_res_eeg = out_of_sample_prediction(trn_data = trn_eeg_data,
val_data = val_eeg_data,
dv_var="fr_tot",
pred_var="value")
print(sprintf("EEG R-sq = %f",oos_res_eeg$rsq))
## [1] "EEG R-sq = 0.969967"
trn_lfp_data = camk_pow_data %>% filter(dset=="training") %>% filter(measure=="LFP")
val_lfp_data = camk_pow_data %>% filter(dset=="validation") %>% filter(measure=="LFP")
oos_res_lfp = out_of_sample_prediction(trn_data = trn_lfp_data,
val_data = val_lfp_data,
dv_var="fr_tot",
pred_var="value")
print(sprintf("LFP R-sq = %f",oos_res_lfp$rsq))
## [1] "LFP R-sq = 0.972305"
# make plot of actual vs predicted
data2plot = rbind(oos_res_eeg$val_data, oos_res_lfp$val_data)
vlineval = NULL
yLabel = "Predicted Firing Rate (Hz)"
xLabel = "Actual Firing Rate (Hz)"
cols2use = dtype_cols2use
plot_diag_line = TRUE
pow_camk_El_oos_plot = make_scatterplot(data2plot = data2plot,
x_var = "fr_tot",
y_var = "pred",
strat_var = "measure",
vlineval = vlineval,
xlabel2use = xLabel,
ylabel2use = yLabel,
title2use = title2use,
cols2use = cols2use,
plot_diag_line = plot_diag_line)
pow_camk_El_oos_plot
Analysis of hM4Di
title2use = "hM4Di"
# ==============================================================================
# H
# Scatterplot against El
data2plot = hm4di_h_data %>% filter(dset=="training")
vlineval = -70
yLabel = "H"
xLabel = "Resting Potential (mV)"
cols2use = dtype_cols2use
plot_fit_line = TRUE
H_hm4di_El_plot = make_scatterplot(data2plot = data2plot,
x_var = "param",
y_var = "value",
strat_var = "measure",
vlineval = vlineval,
xlabel2use = xLabel,
ylabel2use = yLabel,
title2use = title2use,
cols2use = cols2use,
plot_fit_line = plot_fit_line)
H_hm4di_El_plot
# correlation with El in LFP
data2use = data2plot %>% filter(measure=="LFP")
res = cor.test(data2use[,"param"], data2use[,"value"])
res
##
## Pearson's product-moment correlation
##
## data: data2use[, "param"] and data2use[, "value"]
## t = -10.598, df = 9, p-value = 2.203e-06
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## -0.9904076 -0.8569264
## sample estimates:
## cor
## -0.9621892
res$p.value
## [1] 2.202948e-06
# correlation with El in EEG
data2use = data2plot %>% filter(measure=="EEG")
res = cor.test(data2use[,"param"], data2use[,"value"])
res
##
## Pearson's product-moment correlation
##
## data: data2use[, "param"] and data2use[, "value"]
## t = -25.916, df = 9, p-value = 9.144e-10
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## -0.9983368 -0.9737372
## sample estimates:
## cor
## -0.9933664
res$p.value
## [1] 9.144413e-10
# ==============================================================================
# out of sample prediction
trn_eeg_data = hm4di_h_data %>% filter(dset=="training") %>% filter(measure=="EEG")
val_eeg_data = hm4di_h_data %>% filter(dset=="validation") %>% filter(measure=="EEG")
oos_res_eeg = out_of_sample_prediction(trn_data = trn_eeg_data,
val_data = val_eeg_data,
dv_var="fr_tot",
pred_var="value")
print(sprintf("EEG R-sq = %f",oos_res_eeg$rsq))
## [1] "EEG R-sq = 0.919404"
trn_lfp_data = hm4di_h_data %>% filter(dset=="training") %>% filter(measure=="LFP")
val_lfp_data = hm4di_h_data %>% filter(dset=="validation") %>% filter(measure=="LFP")
oos_res_lfp = out_of_sample_prediction(trn_data = trn_lfp_data,
val_data = val_lfp_data,
dv_var="fr_tot",
pred_var="value")
print(sprintf("LFP R-sq = %f",oos_res_lfp$rsq))
## [1] "LFP R-sq = 0.947733"
# make plot of actual vs predicted
data2plot = rbind(oos_res_eeg$val_data, oos_res_lfp$val_data)
vlineval = NULL
yLabel = "Predicted Firing Rate (Hz)"
xLabel = "Actual Firing Rate (Hz)"
cols2use = dtype_cols2use
plot_diag_line = TRUE
H_hm4di_El_oos_plot = make_scatterplot(data2plot = data2plot,
x_var = "fr_tot",
y_var = "pred",
strat_var = "measure",
vlineval = vlineval,
xlabel2use = xLabel,
ylabel2use = yLabel,
title2use = title2use,
cols2use = cols2use,
plot_diag_line = plot_diag_line)
H_hm4di_El_oos_plot
# ==============================================================================
# 1/f slope
# Scatterplot against El
data2plot = hm4di_slope_data %>% filter(dset=="training")
vlineval = -70
yLabel = "1/f slope"
xLabel = "Resting Potential (mV)"
cols2use = dtype_cols2use
plot_fit_line = TRUE
slope_hm4di_El_plot = make_scatterplot(data2plot = data2plot,
x_var = "param",
y_var = "value",
strat_var = "measure",
vlineval = vlineval,
xlabel2use = xLabel,
ylabel2use = yLabel,
title2use = title2use,
cols2use = cols2use,
plot_fit_line = plot_fit_line)
slope_hm4di_El_plot
# correlation with g in LFP
data2use = data2plot %>% filter(measure=="LFP")
res = cor.test(data2use[,"param"], data2use[,"value"])
res
##
## Pearson's product-moment correlation
##
## data: data2use[, "param"] and data2use[, "value"]
## t = 4.5528, df = 9, p-value = 0.00138
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.4711351 0.9560171
## sample estimates:
## cor
## 0.8350161
res$p.value
## [1] 0.001380469
# correlation with g in EEG
data2use = data2plot %>% filter(measure=="EEG")
res = cor.test(data2use[,"param"], data2use[,"value"])
res
##
## Pearson's product-moment correlation
##
## data: data2use[, "param"] and data2use[, "value"]
## t = 4.3696, df = 9, p-value = 0.001798
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.4442079 0.9529888
## sample estimates:
## cor
## 0.8244045
res$p.value
## [1] 0.001798067
# ==============================================================================
# out of sample prediction
trn_eeg_data = hm4di_slope_data %>% filter(dset=="training") %>% filter(measure=="EEG")
val_eeg_data = hm4di_slope_data %>% filter(dset=="validation") %>% filter(measure=="EEG")
oos_res_eeg = out_of_sample_prediction(trn_data = trn_eeg_data,
val_data = val_eeg_data,
dv_var="fr_tot",
pred_var="value")
print(sprintf("EEG R-sq = %f",oos_res_eeg$rsq))
## [1] "EEG R-sq = 0.658198"
trn_lfp_data = hm4di_slope_data %>% filter(dset=="training") %>% filter(measure=="LFP")
val_lfp_data = hm4di_slope_data %>% filter(dset=="validation") %>% filter(measure=="LFP")
oos_res_lfp = out_of_sample_prediction(trn_data = trn_lfp_data,
val_data = val_lfp_data,
dv_var="fr_tot",
pred_var="value")
print(sprintf("LFP R-sq = %f",oos_res_lfp$rsq))
## [1] "LFP R-sq = 0.676430"
# make plot of actual vs predicted
data2plot = rbind(oos_res_eeg$val_data, oos_res_lfp$val_data)
vlineval = NULL
yLabel = "Predicted Firing Rate (Hz)"
xLabel = "Actual Firing Rate (Hz)"
cols2use = dtype_cols2use
plot_diag_line = TRUE
slope_hm4di_El_oos_plot = make_scatterplot(data2plot = data2plot,
x_var = "fr_tot",
y_var = "pred",
strat_var = "measure",
vlineval = vlineval,
xlabel2use = xLabel,
ylabel2use = yLabel,
title2use = title2use,
cols2use = cols2use,
plot_diag_line = plot_diag_line)
slope_hm4di_El_oos_plot
# ==============================================================================
# Power
# Scatterplot against g
data2plot = hm4di_pow_data %>% filter(dset=="training")
vlineval = -70
yLabel = "Broadband Power (a.u.)"
xLabel = "Resting Potential (mV)"
cols2use = dtype_cols2use
plot_fit_line = TRUE
pow_hm4di_El_plot = make_scatterplot(data2plot = data2plot,
x_var = "param",
y_var = "value",
strat_var = "measure",
vlineval = vlineval,
xlabel2use = xLabel,
ylabel2use = yLabel,
title2use = title2use,
cols2use = cols2use,
plot_fit_line = plot_fit_line)
pow_hm4di_El_plot
# correlation with g in LFP
data2use = data2plot %>% filter(measure=="LFP")
res = cor.test(data2use[,"param"], data2use[,"value"])
res
##
## Pearson's product-moment correlation
##
## data: data2use[, "param"] and data2use[, "value"]
## t = 6.5724, df = 9, p-value = 0.0001025
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.6820298 0.9766279
## sample estimates:
## cor
## 0.9097119
res$p.value
## [1] 0.0001024801
# correlation with g in EEG
data2use = data2plot %>% filter(measure=="EEG")
res = cor.test(data2use[,"param"], data2use[,"value"])
res
##
## Pearson's product-moment correlation
##
## data: data2use[, "param"] and data2use[, "value"]
## t = 5.1283, df = 9, p-value = 0.000621
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.5459841 0.9639247
## sample estimates:
## cor
## 0.8631563
res$p.value
## [1] 0.0006209952
# ==============================================================================
# out of sample prediction
trn_eeg_data = hm4di_pow_data %>% filter(dset=="training") %>% filter(measure=="EEG")
val_eeg_data = hm4di_pow_data %>% filter(dset=="validation") %>% filter(measure=="EEG")
oos_res_eeg = out_of_sample_prediction(trn_data = trn_eeg_data,
val_data = val_eeg_data,
dv_var="fr_tot",
pred_var="value")
print(sprintf("EEG R-sq = %f",oos_res_eeg$rsq))
## [1] "EEG R-sq = 0.940921"
trn_lfp_data = hm4di_pow_data %>% filter(dset=="training") %>% filter(measure=="LFP")
val_lfp_data = hm4di_pow_data %>% filter(dset=="validation") %>% filter(measure=="LFP")
oos_res_lfp = out_of_sample_prediction(trn_data = trn_lfp_data,
val_data = val_lfp_data,
dv_var="fr_tot",
pred_var="value")
print(sprintf("LFP R-sq = %f",oos_res_lfp$rsq))
## [1] "LFP R-sq = 0.964205"
# make plot of actual vs predicted
data2plot = rbind(oos_res_eeg$val_data, oos_res_lfp$val_data)
vlineval = NULL
yLabel = "Predicted Firing Rate (Hz)"
xLabel = "Actual Firing Rate (Hz)"
cols2use = dtype_cols2use
plot_diag_line = TRUE
pow_hm4di_El_oos_plot = make_scatterplot(data2plot = data2plot,
x_var = "fr_tot",
y_var = "pred",
strat_var = "measure",
vlineval = vlineval,
xlabel2use = xLabel,
ylabel2use = yLabel,
title2use = title2use,
cols2use = cols2use,
plot_diag_line = plot_diag_line)
pow_hm4di_El_oos_plot
Plots
fig1 = (FR_g_plot + FR_camk_El_plot + FR_hm4di_El_plot) /
(H_g_plot + H_camk_El_plot + H_hm4di_El_plot) /
(H_g_oos_plot + H_camk_El_oos_plot + H_hm4di_El_oos_plot)
ggsave(filename = file.path(plotdir, "Figure01.pdf"), plot = fig1, width = 12, height = 12)
fig1
supp_fig1a = (slope_g_plot + slope_camk_El_plot + slope_hm4di_El_plot) /
(slope_g_oos_plot + slope_camk_El_oos_plot + slope_hm4di_El_oos_plot)
ggsave(filename = file.path(plotdir, "SuppFig01a.pdf"), plot = supp_fig1a, width = 12, height = 8)
supp_fig1a
supp_fig1b = (pow_g_plot + pow_camk_El_plot + pow_hm4di_El_plot) /
(pow_g_oos_plot + pow_camk_El_oos_plot + pow_hm4di_El_oos_plot)
ggsave(filename = file.path(plotdir, "SuppFig01b.pdf"), plot = supp_fig1b, width = 12, height = 8)
supp_fig1b
